clear;close all;filebrowser;workspace;warning('off','all');
h(1) = figure('units','normalized','WindowStyle','docked','outerposition',[0 0 1 1]);
% h(4) = figure('units','normalized','WindowStyle','docked','outerposition',[0 0 1 1]);
% h(2) = figure('units','normalized','WindowStyle','docked','outerposition',[0 0 1 1]);
% h(3) = figure('units','normalized','WindowStyle','docked','outerposition',[0 0 1 1]);
hc = 1.23984193;    %[eV*um]
eV2cm1 = 8065.54429;    %[cm-1/eV]
VR = 1000;          % [V] repeller voltage
AVMI = 46.7;        %[mm/sqrt(meV/V)]
LabelFontSize = 14;

%%--------ODT effect with 356 nm data (5/24/2019)---------------
figure(h(1));

tODToff_unsort    = [0.5   10.5   20.5   2.5   70.5   35.5   50.5]; %[nm]
Ncycle_unsort     = [202   200    200    200   200    200    200];
UVpower_unsort    = 57.*ones(size(tODToff_unsort));   %[mW]
UVwx_unsort       = 345.*ones(size(tODToff_unsort));  %[um]
UVwy_unsort       = 348.*ones(size(tODToff_unsort));  %[um]
UVpower_unsort    = UVpower_unsort./UVwx_unsort./UVwy_unsort*350^2;


NK_unsort         = [1.2   1.2   1.2    1.3   1.2     1.2    1.1]*1e3;
NKbkg_unsort      = [21    35    24     24    23      33     22];

NRb_unsort        = [3.5   3.9   3.9    4.0   3.5     4.5    3.4]*1e3;
NRbbkg_unsort     = [18    23    22     15    12      20     13];

NKRb_unsort       = [3.6   3.2   3.3    3.7   3.8     3.2    3.3]*1e3;
NKRbbkg_unsort    = [19    10    10     21    18      24     11];

NK2_unsort        = [24    26    20     23    20      29     28];
NK2bkg_unsort     = [20    15    13     23    23      12     18];

NRb2_unsort       = [42    21    32     32    23      26     30];
NRb2bkg_unsort    = [14    12    17     13    14      9      19];

NK2Rb_unsort      = [15    7     9      16    12      11     16];
NK2Rbbkg_unsort   = [14    26    19     12    15      15     15];

NKRb2_unsort      = [17    65    128    37    529     175    320];
NKRb2bkg_unsort   = [16    14    13     7     19      10     14];

% NK2Rb2_unsort     = [85    124   108    99    81      120    103];
% NK2Rb2bkg_unsort  = [9     7     8      9     16      18     11];

NK2Rb2_unsort     = [80    110   99    85     69      107    92];
NK2Rb2bkg_unsort  = [9     7     8      9     16      18     11];


% NKRb2norm_unsort  = [];
% NK2Rb2norm_unsort = [];

[tODToff, tODToff_order] = sort(tODToff_unsort);
Ncycle = Ncycle_unsort(tODToff_order);
UVpower = UVpower_unsort(tODToff_order);

NK   =   NK_unsort(tODToff_order);
NKbkg = NKbkg_unsort(tODToff_order);
NRb  =   NRb_unsort(tODToff_order);
NRbbkg = NRbbkg_unsort(tODToff_order);

NKRb =   NKRb_unsort(tODToff_order);
NKRbbkg =   NKRbbkg_unsort(tODToff_order);
NK2  =   NK2_unsort(tODToff_order);
NK2bkg = NK2bkg_unsort(tODToff_order);
NRb2 =   NRb2_unsort(tODToff_order);
NRb2bkg = NRb2bkg_unsort(tODToff_order);
NK2Rb = NK2Rb_unsort(tODToff_order);
NK2Rbbkg = NK2Rbbkg_unsort(tODToff_order);
NKRb2 = NKRb2_unsort(tODToff_order);
NKRb2bkg = NKRb2bkg_unsort(tODToff_order);
NK2Rb2 = NK2Rb2_unsort(tODToff_order);
NK2Rb2bkg = NK2Rb2bkg_unsort(tODToff_order);
% NKRb2norm = NKRb2norm_unsort(tODToff_order);
% NK2Rb2norm = NK2Rb2norm_unsort(tODToff_order);

markersize = 5;
subplot(2, 2, 1)
err = sqrt(NK2Rb2+(0.1*NK2Rb2).^2);
errorbar(tODToff, NK2Rb2, err, 'ob', 'MarkerSize', markersize,...
    'MarkerFaceColor', 'b', ...
    'color','k', 'CapSize', 1, 'LineWidth', 0.5);
%----------linear fit-----------
linearFit = @(b, x) b(1).*x + b(2);
startPoints = [0 100];
W = err'.^(-2);
[b, residual, jacobian,~,~,~] = nlinfit(tODToff', NK2Rb2', linearFit, startPoints, 'Weights', W)
%%calculate 95% confidence interval of the fitted parameters
alpha = 0.05;       %95% level
ci = nlparci(b, residual, 'Jacobian', jacobian, 'alpha', alpha);
b = b';
b_Err = b - ci(:, 1);
%%calculate SE
t = tinv(1-alpha/2,length(tODToff)-length(b));
nlinfit_se = (ci(:,2)-ci(:,1)) ./ (2*t);  % Standard Error

disp('--------fit result-------');
disp(['Slope = ', num2str(b(1)), ' +/- ', num2str(nlinfit_se(1)),...
    ', 95% level: +/-', num2str(b_Err(1)), '.']);
disp(['y0 = ', num2str(b(2)), ' +/- ', num2str(nlinfit_se(2)),...
    ', 95% level: +/-', num2str(b_Err(2)), '.']);
hold on
xplot = 0:75;
yplot = linearFit(b,xplot);
h2 = plot(xplot, yplot);
set(h2, 'LineWidth', 4, 'Color', 'b');
h2.Color(4)=0.5;  % 70% transparent
legend off
hold off
% errorbar(EUV, RVMIKRb2.^2, 2.*RVMIKRb2.*RVMIKRb2_err, ...
%     'ok', 'MarkerSize',7,'MarkerEdgeColor',[0.5,0.5,0.5], ...
%     'MarkerFaceColor','m', 'color','k', 'CapSize', 1, 'LineWidth', 0.5);

% hold on
% plot(tODToff, (NKRb2-NKRb2bkg)./Ncycle*200./UVpower*245, '-om', 'MarkerSize', markersize,'MarkerFaceColor', 'm');
% plot(tODToff, (NK2-NK2bkg)./Ncycle*200./UVpower*245, '-or');
% plot(tODToff, (NRb2-NRb2bkg)./Ncycle*200./UVpower*245, '-oc', 'MarkerSize', markersize);
% plot(tODToff, (NK2Rb-NK2Rbbkg)./Ncycle*200./UVpower*245, '-og', 'MarkerSize', markersize);
xlabel('t_{off} (\mus)');
xlim([-3, 75]);
ylabel('Ion (count)');
title('ODT effect with 356 nm UV ionization');
ylim([0, 150]);
hold off
% legend('(N_{K_2Rb_2}-N_{K_2Rb_2bkg})',...
%     '(N_{KRb_2}-N_{KRb_2bkg})',...
%     '(N_{K_2}-N_{K_2bkg})',...
%     '(N_{Rb_2}-N_{Rb_2bkg})',...
%     '(N_{K_2Rb}-N_{K_2Rb bkg})',...
%     'Location', 'northwest');
subplot(2,2, 3)
plot(tODToff, (NK-NKbkg)./Ncycle*200./UVpower*245, '-ok', 'MarkerSize', markersize,'MarkerFaceColor', 'b');
hold on;
plot(tODToff, (NRb-NRbbkg)./Ncycle*200./UVpower*245, '-ok', 'MarkerSize', markersize,'MarkerFaceColor', 'r');
plot(tODToff, (NKRb-NKRbbkg)./Ncycle*200./UVpower*245, '-ok', 'MarkerSize', markersize,'MarkerFaceColor', 'k');
hold off;
xlabel('ODT off time before ionization (\mus)');
ylabel('Ion (count)');
ylim([0 7]*1e4);
legend('(N_{K}-N_{K bkg})',...
    '(N_{Rb}-N_{Rb bkg})',...
    '(N_{KRb}-N_{KRb bkg})',...
    'Location', 'northeast');
disp('=======================');



%%--------ODT effect with 340 nm data (5/24/2019)---------------
figure(h(1));

tODToff_unsort    = [50.5  0.5   10.5   20.5    35.5   70.5]; %[nm]
Ncycle_unsort     = [200   200   200    200     200    200];
UVpower_unsort    = 204.*ones(size(tODToff_unsort));   %[mW]
UVwx_unsort       = 385.*ones(size(tODToff_unsort));  %[um]
UVwy_unsort       = 362.*ones(size(tODToff_unsort));  %[um]
UVpower_unsort    = UVpower_unsort./UVwx_unsort./UVwy_unsort*350^2;


NK_unsort         = [38.9  38.7  41.1   43.7   44.7   45.7]*1e3;
NKbkg_unsort      = [52    43    38     41     33     44 ];

NRb_unsort        = [8.9   7.1   8.3    8.9    9.6    9.5]*1e3;
NRbbkg_unsort     = [19    18    19     14     22     22];

NKRb_unsort       = [8.8   8.4   8.8    8.9    9.1    9.2]*1e3;
NKRbbkg_unsort    = [18    16    16     26     17     12];

NK2_unsort        = [18    19    23     26     18     14];
NK2bkg_unsort     = [20    23    16     18     32     15];

NRb2_unsort       = [52    44    40     46     47     40];
NRb2bkg_unsort    = [11    16    11     12     12     9];

NK2Rb_unsort      = [14    20    15     15     16     14];
NK2Rbbkg_unsort   = [18    11    14     13     11     20];

NKRb2_unsort      = [214   28    106    148    188    221];
NKRb2bkg_unsort   = [9     9     16     16     10     12];

NK2Rb2_unsort     = [34    17    29     34     31     26];
NK2Rb2bkg_unsort  = [23    19    19     18     19     22];

% NKRb2norm_unsort  = [];
% NK2Rb2norm_unsort = [];

[tODToff, tODToff_order] = sort(tODToff_unsort);
Ncycle = Ncycle_unsort(tODToff_order);
UVpower = UVpower_unsort(tODToff_order);

NK   =   NK_unsort(tODToff_order);
NKbkg = NKbkg_unsort(tODToff_order);
NRb  =   NRb_unsort(tODToff_order);
NRbbkg = NRbbkg_unsort(tODToff_order);

NKRb =   NKRb_unsort(tODToff_order);
NKRbbkg =   NKRbbkg_unsort(tODToff_order);
NK2  =   NK2_unsort(tODToff_order);
NK2bkg = NK2bkg_unsort(tODToff_order);
NRb2 =   NRb2_unsort(tODToff_order);
NRb2bkg = NRb2bkg_unsort(tODToff_order);
NK2Rb = NK2Rb_unsort(tODToff_order);
NK2Rbbkg = NK2Rbbkg_unsort(tODToff_order);
NKRb2 = NKRb2_unsort(tODToff_order);
NKRb2bkg = NKRb2bkg_unsort(tODToff_order);
NK2Rb2 = NK2Rb2_unsort(tODToff_order);
NK2Rb2bkg = NK2Rb2bkg_unsort(tODToff_order);
% NKRb2norm = NKRb2norm_unsort(tODToff_order);
% NK2Rb2norm = NK2Rb2norm_unsort(tODToff_order);

markersize = 5;
subplot(2,2,2)
plot(tODToff, (NK2Rb2-NK2Rb2bkg)./Ncycle*200./UVpower*245, '-ob', 'MarkerSize', markersize,'MarkerFaceColor', 'b');
hold on
plot(tODToff, (NKRb2-NKRb2bkg)./Ncycle*200./UVpower*245, '-om', 'MarkerSize', markersize,'MarkerFaceColor', 'm');
% plot(tODToff, (NK2-NK2bkg)./Ncycle*200./UVpower*245, '-or');
% plot(tODToff, (NRb2-NRb2bkg)./Ncycle*200./UVpower*245, '-oc', 'MarkerSize', markersize);
% plot(tODToff, (NK2Rb-NK2Rbbkg)./Ncycle*200./UVpower*245, '-og', 'MarkerSize', markersize);
xlabel('ODT off time before ionization (\mus)');
ylabel('Ion (count)');
title('ODT effect with 340 nm UV ionization');
ylim([-200, 2300]);
hold off
legend('(N_{K_2Rb_2}-N_{K_2Rb_2bkg})',...
    '(N_{KRb_2}-N_{KRb_2bkg})',...
    '(N_{K_2}-N_{K_2bkg})',...
    '(N_{Rb_2}-N_{Rb_2bkg})',...
    '(N_{K_2Rb}-N_{K_2Rb bkg})',...
    'Location', 'northwest');
subplot(2,2,4)
plot(tODToff, (NK-NKbkg)./Ncycle*200./UVpower*245, '-ok', 'MarkerSize', markersize,'MarkerFaceColor', 'b');
hold on;
plot(tODToff, (NRb-NRbbkg)./Ncycle*200./UVpower*245, '-ok', 'MarkerSize', markersize,'MarkerFaceColor', 'r');
plot(tODToff, (NKRb-NKRbbkg)./Ncycle*200./UVpower*245, '-ok', 'MarkerSize', markersize,'MarkerFaceColor', 'k');
hold off;
xlabel('ODT off time before ionization (\mus)');
ylabel('Ion (count)');
ylim([0 7]*1e4);
% legend('(N_{K}-N_{K bkg})',...
%     '(N_{Rb}-N_{Rb bkg})',...
%     '(N_{KRb}-N_{KRb bkg})',...
%     'Location', 'northeastoutside');
disp('=======================');

%%%-------save figures--------------
savefig(h, ['plots/summaryplot_ODTeffect_2.fig']);
saveas(h(1), ['plots/summaryplot_ODTeffect_2.png']);
% saveas(h(3), ['plots/summaryplot_K2Rb2_spectroscopy.png']);
